#######################################################################
## author:    Michael L. Wicki
## contact:   wimi@ethz.ch, ETH Zurich
## file name: conjoint_EVZH_english.R
## Context:   ISTP mobility lab - EVZH
## started:   2021-03-11
## Summary:   Analysis Survey Conjoint EVZH
#######################################################################

rm(list = ls())
getwd()
setwd("Z:/02_People/Wicki Michael/EVZH")
set.seed(42)

library(cjoint)
library(Matrix)
library(ggthemes)
library(ggplot2)
library(data.table)
library(dplyr)
library(nFactors)
library(lme4)
library(effects)
library(ggthemes)
library(tidyr)
library(forcats)
library(texreg)
library(influence.ME)

########### load conjoint data ###############
load("Z:/02_People/Wicki Michael/EVZH/evzh_conjoint.RData")


######################### create factors with baselines #########################

df_conj$Kaufpreis <- factor(df_conj$Kaufpreis, levels=c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000", "CHF 100'000"))
df_conj$Energiekosten <- factor(df_conj$Energiekosten, levels=c("CHF 2",  "CHF 5", "CHF 10","CHF 20", "CHF 30"))
df_conj$Tageskarte <- factor(df_conj$Tageskarte, levels=c("CHF 0",  "CHF 25",  "CHF 75", "CHF 150"))
df_conj$Verbrennertausch <- factor(df_conj$Verbrennertausch, levels=c("CHF 0", "CHF 50", "CHF 100", "CHF 150"))
df_conj$Beitrag <- factor(df_conj$Beitrag, levels=c("CHF 0", "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000"))
df_conj$Ladedauer <- factor(df_conj$Ladedauer, levels=c("15 min", "30 min", "45 min", "60 min", "75 min", "90 min"))
df_conj$Unterhaltskosten <- factor(df_conj$Unterhaltskosten, levels=c("CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000"))
df_conj$Garantie <- factor(df_conj$Garantie, levels=c("2 Jahre",  "4 Jahre",  "6 Jahre",  "8 Jahre", "10 Jahre"))
df_conj$Reichweite <- factor(df_conj$Reichweite, levels=c("100 km", "300 km",  "500 km", "1000 km"))
df_conj$Ladegarantie <- factor(df_conj$Ladegarantie, levels=c("0,5 km", "1 km", "2 km", "10 km"))

######################## Recode ############################
baselines_base <- list()
baselines_base$Kaufpreis <- "CHF 20'000"
baselines_base$Energiekosten <- "CHF 2"
baselines_base$Tageskarte <- "CHF 0"
baselines_base$Verbrennertausch <- "CHF 0"
baselines_base$Beitrag <- "CHF 0"
baselines_base$Ladedauer <- "15 min"
baselines_base$Unterhaltskosten <- "CHF 500"
baselines_base$Garantie <- "2 Jahre"
baselines_base$Reichweite <- "100 km"
baselines_base$Ladegarantie <- "0,5 km"

attribute_list <- list()
attribute_list$Kaufpreis <- c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000", "CHF 100'000")
attribute_list$Energiekosten <- c("CHF 2",  "CHF 5", "CHF 10","CHF 20", "CHF 30")
attribute_list$Tageskarte <- c("CHF 0",  "CHF 25",  "CHF 75", "CHF 150")
attribute_list$Verbrennertausch <- c("CHF 0", "CHF 50", "CHF 100", "CHF 150")
attribute_list$Beitrag <- c("CHF 0", "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000")
attribute_list$Ladedauer <- c("15 min", "30 min", "45 min", "60 min", "75 min", "90 min")
attribute_list$Unterhaltskosten <- c("CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000")
attribute_list$Garantie <- c("2 Jahre",  "4 Jahre",  "6 Jahre",  "8 Jahre", "10 Jahre")
attribute_list$Reichweite <- c("100 km", "300 km",  "500 km", "1000 km")
attribute_list$Ladegarantie <- c("0,5 km", "1 km", "2 km", "10 km")


constraint_list <- list()

design <- makeDesign(type='constraints', attribute.levels=attribute_list,
                     constraints=constraint_list)


######################### Conjoint Choice All ######################### 
m_base_choice <- lmer(choice ~ Kaufpreis +  Energiekosten +  Tageskarte +  Verbrennertausch + Beitrag +  
                        Ladedauer + Unterhaltskosten + Garantie + Reichweite + Ladegarantie + (1|ResponseId), data = df_conj)

m_base_rate <- lmer(buy ~ Kaufpreis +  Energiekosten +  Tageskarte +  Verbrennertausch + Beitrag +  
                      Ladedauer + Unterhaltskosten + Garantie + Reichweite + Ladegarantie + (1|ResponseId), data = df_conj)


######################### Choice ######################### 
c_allC <- cjoint::amce(choice ~ Kaufpreis +  Energiekosten +  Tageskarte +  Verbrennertausch + Beitrag +  
                         Ladedauer + Unterhaltskosten + Garantie + Reichweite + Ladegarantie, data = df_conj, 
                       cluster=TRUE, respondent.id="ResponseId",
                       design=design, baselines = baselines_base)
summary(c_allC)
plot(c_allC)



c_attribute1 <- c(0, c_allC$estimates$Kaufpreis["AMCE","KaufpreisCHF40000"], c_allC$estimates$Kaufpreis["AMCE","KaufpreisCHF60000"], c_allC$estimates$Kaufpreis["AMCE","KaufpreisCHF80000"], c_allC$estimates$Kaufpreis["AMCE","KaufpreisCHF100000"])
s_attribute1 <- c(0, c_allC$estimates$Kaufpreis["Std. Error","KaufpreisCHF40000"], c_allC$estimates$Kaufpreis["Std. Error","KaufpreisCHF60000"], c_allC$estimates$Kaufpreis["Std. Error","KaufpreisCHF80000"], c_allC$estimates$Kaufpreis["Std. Error","KaufpreisCHF100000"])

c_attribute2 <- c(0, c_allC$estimates$Energiekosten["AMCE","EnergiekostenCHF5"],c_allC$estimates$Energiekosten["AMCE","EnergiekostenCHF10"],c_allC$estimates$Energiekosten["AMCE","EnergiekostenCHF20"],c_allC$estimates$Energiekosten["AMCE","EnergiekostenCHF30"])
s_attribute2 <- c(0, c_allC$estimates$Energiekosten["Std. Error","EnergiekostenCHF5"],c_allC$estimates$Energiekosten["Std. Error","EnergiekostenCHF10"],c_allC$estimates$Energiekosten["Std. Error","EnergiekostenCHF20"],c_allC$estimates$Energiekosten["Std. Error","EnergiekostenCHF30"])

c_attribute3 <- c(0, c_allC$estimates$Tageskarte["AMCE","TageskarteCHF25"], c_allC$estimates$Tageskarte["AMCE","TageskarteCHF75"], c_allC$estimates$Tageskarte["AMCE","TageskarteCHF150"])
s_attribute3 <- c(0, c_allC$estimates$Tageskarte["Std. Error","TageskarteCHF25"], c_allC$estimates$Tageskarte["Std. Error","TageskarteCHF75"], c_allC$estimates$Tageskarte["Std. Error","TageskarteCHF150"])

c_attribute4 <- c(0, c_allC$estimates$Verbrennertausch["AMCE","VerbrennertauschCHF50"], c_allC$estimates$Verbrennertausch["AMCE","VerbrennertauschCHF100"], c_allC$estimates$Verbrennertausch["AMCE","VerbrennertauschCHF150"])
s_attribute4 <- c(0, c_allC$estimates$Verbrennertausch["Std. Error","VerbrennertauschCHF50"], c_allC$estimates$Verbrennertausch["Std. Error","VerbrennertauschCHF100"], c_allC$estimates$Verbrennertausch["Std. Error","VerbrennertauschCHF150"])

c_attribute5 <- c(0, c_allC$estimates$Beitrag["AMCE","BeitragCHF500"], c_allC$estimates$Beitrag["AMCE","BeitragCHF2000"], c_allC$estimates$Beitrag["AMCE","BeitragCHF4000"], c_allC$estimates$Beitrag["AMCE","BeitragCHF6000"])
s_attribute5 <- c(0, c_allC$estimates$Beitrag["Std. Error","BeitragCHF500"], c_allC$estimates$Beitrag["Std. Error","BeitragCHF2000"], c_allC$estimates$Beitrag["Std. Error","BeitragCHF4000"], c_allC$estimates$Beitrag["Std. Error","BeitragCHF6000"])

c_attribute6 <- c(0, c_allC$estimates$Ladedauer["AMCE","Ladedauer30min"], c_allC$estimates$Ladedauer["AMCE","Ladedauer45min"], c_allC$estimates$Ladedauer["AMCE","Ladedauer60min"], c_allC$estimates$Ladedauer["AMCE","Ladedauer75min"], c_allC$estimates$Ladedauer["AMCE","Ladedauer90min"])
s_attribute6 <- c(0, c_allC$estimates$Ladedauer["Std. Error","Ladedauer30min"], c_allC$estimates$Ladedauer["Std. Error","Ladedauer45min"], c_allC$estimates$Ladedauer["Std. Error","Ladedauer60min"], c_allC$estimates$Ladedauer["Std. Error","Ladedauer75min"], c_allC$estimates$Ladedauer["Std. Error","Ladedauer90min"])

c_attribute7 <- c(0, c_allC$estimates$Unterhaltskosten["AMCE","UnterhaltskostenCHF2000"], c_allC$estimates$Unterhaltskosten["AMCE","UnterhaltskostenCHF4000"], c_allC$estimates$Unterhaltskosten["AMCE","UnterhaltskostenCHF6000"])
s_attribute7 <- c(0, c_allC$estimates$Unterhaltskosten["Std. Error","UnterhaltskostenCHF2000"], c_allC$estimates$Unterhaltskosten["Std. Error","UnterhaltskostenCHF4000"], c_allC$estimates$Unterhaltskosten["Std. Error","UnterhaltskostenCHF6000"])

c_attribute8 <- c(0, c_allC$estimates$Garantie["AMCE","Garantie4Jahre"], c_allC$estimates$Garantie["AMCE","Garantie6Jahre"], c_allC$estimates$Garantie["AMCE","Garantie8Jahre"], c_allC$estimates$Garantie["AMCE","Garantie10Jahre"])
s_attribute8 <- c(0, c_allC$estimates$Garantie["Std. Error","Garantie4Jahre"], c_allC$estimates$Garantie["Std. Error","Garantie6Jahre"], c_allC$estimates$Garantie["Std. Error","Garantie8Jahre"], c_allC$estimates$Garantie["Std. Error","Garantie10Jahre"])

c_attribute9 <- c(0, c_allC$estimates$Reichweite["AMCE","Reichweite300km"], c_allC$estimates$Reichweite["AMCE","Reichweite500km"], c_allC$estimates$Reichweite["AMCE","Reichweite1000km"])
s_attribute9 <- c(0, c_allC$estimates$Reichweite["Std. Error","Reichweite300km"], c_allC$estimates$Reichweite["Std. Error","Reichweite500km"], c_allC$estimates$Reichweite["Std. Error","Reichweite1000km"])

c_attribute10 <- c(0, c_allC$estimates$Ladegarantie["AMCE","Ladegarantie1km"], c_allC$estimates$Ladegarantie["AMCE","Ladegarantie2km"], c_allC$estimates$Ladegarantie["AMCE","Ladegarantie10km"])
s_attribute10 <- c(0, c_allC$estimates$Ladegarantie["Std. Error","Ladegarantie1km"], c_allC$estimates$Ladegarantie["Std. Error","Ladegarantie2km"], c_allC$estimates$Ladegarantie["Std. Error","Ladegarantie10km"])


c <- c(c_attribute1, c_attribute2, c_attribute3, c_attribute4, c_attribute5, c_attribute6, c_attribute7, c_attribute8, c_attribute9, c_attribute10)
s <- c(s_attribute1, s_attribute2, s_attribute3, s_attribute4, s_attribute5, s_attribute6, s_attribute7, s_attribute8, s_attribute9, s_attribute10)

df_c_allC <- data.frame(c,s)

c_allC$attributes$Tageskarte <- ifelse(c_allC$attributes$Tageskarte == "CHF0", "1CHF0", ifelse(
  c_allC$attributes$Tageskarte == "CHF150", "1CHF150", c_allC$attributes$Tageskarte))
c_allC$attributes$Beitrag <- ifelse(c_allC$attributes$Beitrag == "CHF0", "2CHF0", ifelse(
  c_allC$attributes$Beitrag == "CHF500", "2CHF500", ifelse(
    c_allC$attributes$Beitrag == "CHF2000", "2CHF2000", ifelse(
      c_allC$attributes$Beitrag == "CHF4000", "2CHF4000", ifelse(
        c_allC$attributes$Beitrag == "CHF6000", "2CHF6000", c_allC$attributes$Beitrag)))))


df_c_allC$att <- c(c_allC$attributes$Kaufpreis[1], c_allC$attributes$Kaufpreis[2:5], 
                   c_allC$attributes$Energiekosten[1], c_allC$attributes$Energiekosten[2:5],
                   c_allC$attributes$Tageskarte[1], c_allC$attributes$Tageskarte[2:4], 
                   c_allC$attributes$Verbrennertausch[1], c_allC$attributes$Verbrennertausch[2:4],
                   c_allC$attributes$Beitrag[1], c_allC$attributes$Beitrag[2:5], 
                   c_allC$attributes$Ladedauer[1], c_allC$attributes$Ladedauer[2:6],
                   c_allC$attributes$Unterhaltskosten[1], c_allC$attributes$Unterhaltskosten[2:4],
                   c_allC$attributes$Garantie[1], c_allC$attributes$Garantie[2:5],
                   c_allC$attributes$Reichweite[1], c_allC$attributes$Reichweite[2:4],
                   c_allC$attributes$Ladegarantie[1], c_allC$attributes$Ladegarantie[2:4])

df_c_allC$policy <- ifelse(df_c_allC$att == "CHF20000" |  df_c_allC$att == "CHF40000" |  df_c_allC$att == "CHF60000" |  df_c_allC$att == "CHF80000" |  df_c_allC$att == "CHF100000", 1, ifelse(
  df_c_allC$att == "CHF2" | df_c_allC$att == "CHF5" | df_c_allC$att == "CHF10" | df_c_allC$att == "CHF20" | df_c_allC$att == "CHF30", 2, ifelse(
    df_c_allC$att == "1CHF0" | df_c_allC$att == "CHF25" |  df_c_allC$att == "CHF75" | df_c_allC$att == "1CHF150", 3, ifelse(
      df_c_allC$att == "CHF0" | df_c_allC$att == "CHF50" | df_c_allC$att == "CHF100" | df_c_allC$att == "CHF150", 4, ifelse(
        df_c_allC$att == "2CHF0" | df_c_allC$att == "2CHF500" | df_c_allC$att == "2CHF2000" | df_c_allC$att == "2CHF4000" | df_c_allC$att == "2CHF6000", 5, ifelse(
          df_c_allC$att == "15min" | df_c_allC$att == "30min" | df_c_allC$att == "45min" | df_c_allC$att == "60min" | df_c_allC$att == "75min" | df_c_allC$att == "90min", 6, ifelse(
            df_c_allC$att == "CHF500" | df_c_allC$att == "CHF2000" | df_c_allC$att == "CHF4000" | df_c_allC$att == "CHF6000", 7, ifelse(
            df_c_allC$att == "2Jahre" | df_c_allC$att == "4Jahre" | df_c_allC$att == "6Jahre" | df_c_allC$att == "8Jahre" | df_c_allC$att == "10Jahre", 8, ifelse(
              df_c_allC$att == "100km" | df_c_allC$att == "300km" | df_c_allC$att == "500km" | df_c_allC$att == "1000km", 9, ifelse(
                  df_c_allC$att == "05km" | df_c_allC$att == "1km" | df_c_allC$att == "2km" | df_c_allC$att == "10km", 10, 11))))))))))

df_c_allC$policy = factor(as.factor(df_c_allC$policy), levels=c(1:11), labels=c("Purchase\nprice","Energy\ncosts","PT day\npass", "Car\nexchange", 
                                                                                "Government\nsubsidy", "Charging\ntime", "Maintenance\ncosts", "Warranty\nperiod",
                                                                                "Driving\nrange", "Charging\nguarantee", "Other"))


df_c_allC$att <- factor(df_c_allC$att, levels = c(rev(c("CHF20000","CHF40000","CHF60000","CHF80000","CHF100000",
                                                        "CHF2","CHF5","CHF10", "CHF20", "CHF30", 
                                                        "1CHF0","CHF25", "CHF75", "1CHF150", 
                                                        "CHF0","CHF50", "CHF100","CHF150",
                                                        "2CHF0","2CHF500","2CHF2000", "2CHF4000", "2CHF6000", 
                                                        "15min", "30min", "45min", "60min", "75min","90min", 
                                                        "CHF500","CHF2000", "CHF4000", "CHF6000",
                                                        "2Jahre","4Jahre","6Jahre","8Jahre","10Jahre",
                                                        "100km", "300km", "500km", "1000km",
                                                        "05km","1km", "2km", "10km"))),
                        labels = c(rev(c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000", "CHF 100'000",
                                         "CHF 2",  "CHF 5", "CHF 10","CHF 20", "CHF 30",
                                         "CHF 0",  "CHF 25",  "CHF 75", "CHF 150",
                                         " CHF 0", "CHF 50", "CHF 100", "CHF 150",
                                         "  CHF 0", "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                         "15 min", "30 min", "45 min", "60 min", "75 min", "90 min",
                                         "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                         "2 years",  "4 years",  "6 years",  "8 years", "10 years",
                                         "100 km", "300 km",  "500 km", "1'000 km",
                                         "0,5 km", "1 km", "2 km", "10 km"))))




figall <- ggplot(df_c_allC, aes(x=att, y = c)) +
  geom_pointrange(aes(min = c - 1.95 * s, max = c + 1.95 * s)) +
  theme_bw() + 
  theme(panel.grid.major.x = element_line(size=.3, color="gray80"),
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(size=.3, color="gray80" ),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom",
        strip.background = element_rect(fill = 'white')) +
  coord_flip() +
  ylab("Estimated AMCE") +
  xlab("Attributes") +
  geom_hline(yintercept = 0, lty="dashed") +
  facet_grid(policy ~ ., scales = "free") +
  ggtitle("Conjoint choice")
figall

ggsave(figall, filename = "./FigureChoice.png", width = 10.5, height = 7.5, unit="in", dpi = 600)





######################### Buy ######################### 
b_base <- cjoint::amce(buy ~ Kaufpreis +  Energiekosten +  Tageskarte +  Verbrennertausch + Beitrag +  
                         Ladedauer + Unterhaltskosten + Garantie + Reichweite + Ladegarantie, data = df_conj, 
                       cluster=TRUE, respondent.id="ResponseId",
                       design=design, baselines = baselines_base)
plot(b_base)

c_attribute1 <- c(0, b_base$estimates$Kaufpreis["AMCE","KaufpreisCHF40000"], b_base$estimates$Kaufpreis["AMCE","KaufpreisCHF60000"], b_base$estimates$Kaufpreis["AMCE","KaufpreisCHF80000"], b_base$estimates$Kaufpreis["AMCE","KaufpreisCHF100000"])
s_attribute1 <- c(0, b_base$estimates$Kaufpreis["Std. Error","KaufpreisCHF40000"], b_base$estimates$Kaufpreis["Std. Error","KaufpreisCHF60000"], b_base$estimates$Kaufpreis["Std. Error","KaufpreisCHF80000"], b_base$estimates$Kaufpreis["Std. Error","KaufpreisCHF100000"])

c_attribute2 <- c(0, b_base$estimates$Energiekosten["AMCE","EnergiekostenCHF5"],b_base$estimates$Energiekosten["AMCE","EnergiekostenCHF10"],b_base$estimates$Energiekosten["AMCE","EnergiekostenCHF20"],b_base$estimates$Energiekosten["AMCE","EnergiekostenCHF30"])
s_attribute2 <- c(0, b_base$estimates$Energiekosten["Std. Error","EnergiekostenCHF5"],b_base$estimates$Energiekosten["Std. Error","EnergiekostenCHF10"],b_base$estimates$Energiekosten["Std. Error","EnergiekostenCHF20"],b_base$estimates$Energiekosten["Std. Error","EnergiekostenCHF30"])

c_attribute3 <- c(0, b_base$estimates$Tageskarte["AMCE","TageskarteCHF25"], b_base$estimates$Tageskarte["AMCE","TageskarteCHF75"], b_base$estimates$Tageskarte["AMCE","TageskarteCHF150"])
s_attribute3 <- c(0, b_base$estimates$Tageskarte["Std. Error","TageskarteCHF25"], b_base$estimates$Tageskarte["Std. Error","TageskarteCHF75"], b_base$estimates$Tageskarte["Std. Error","TageskarteCHF150"])

c_attribute4 <- c(0, b_base$estimates$Verbrennertausch["AMCE","VerbrennertauschCHF50"], b_base$estimates$Verbrennertausch["AMCE","VerbrennertauschCHF100"], b_base$estimates$Verbrennertausch["AMCE","VerbrennertauschCHF150"])
s_attribute4 <- c(0, b_base$estimates$Verbrennertausch["Std. Error","VerbrennertauschCHF50"], b_base$estimates$Verbrennertausch["Std. Error","VerbrennertauschCHF100"], b_base$estimates$Verbrennertausch["Std. Error","VerbrennertauschCHF150"])

c_attribute5 <- c(0, b_base$estimates$Beitrag["AMCE","BeitragCHF500"], b_base$estimates$Beitrag["AMCE","BeitragCHF2000"], b_base$estimates$Beitrag["AMCE","BeitragCHF4000"], b_base$estimates$Beitrag["AMCE","BeitragCHF6000"])
s_attribute5 <- c(0, b_base$estimates$Beitrag["Std. Error","BeitragCHF500"], b_base$estimates$Beitrag["Std. Error","BeitragCHF2000"], b_base$estimates$Beitrag["Std. Error","BeitragCHF4000"], b_base$estimates$Beitrag["Std. Error","BeitragCHF6000"])

c_attribute6 <- c(0, b_base$estimates$Ladedauer["AMCE","Ladedauer30min"], b_base$estimates$Ladedauer["AMCE","Ladedauer45min"], b_base$estimates$Ladedauer["AMCE","Ladedauer60min"], b_base$estimates$Ladedauer["AMCE","Ladedauer75min"], b_base$estimates$Ladedauer["AMCE","Ladedauer90min"])
s_attribute6 <- c(0, b_base$estimates$Ladedauer["Std. Error","Ladedauer30min"], b_base$estimates$Ladedauer["Std. Error","Ladedauer45min"], b_base$estimates$Ladedauer["Std. Error","Ladedauer60min"], b_base$estimates$Ladedauer["Std. Error","Ladedauer75min"], b_base$estimates$Ladedauer["Std. Error","Ladedauer90min"])

c_attribute7 <- c(0, b_base$estimates$Unterhaltskosten["AMCE","UnterhaltskostenCHF2000"], b_base$estimates$Unterhaltskosten["AMCE","UnterhaltskostenCHF4000"], b_base$estimates$Unterhaltskosten["AMCE","UnterhaltskostenCHF6000"])
s_attribute7 <- c(0, b_base$estimates$Unterhaltskosten["Std. Error","UnterhaltskostenCHF2000"], b_base$estimates$Unterhaltskosten["Std. Error","UnterhaltskostenCHF4000"], b_base$estimates$Unterhaltskosten["Std. Error","UnterhaltskostenCHF6000"])

c_attribute8 <- c(0, b_base$estimates$Garantie["AMCE","Garantie4Jahre"], b_base$estimates$Garantie["AMCE","Garantie6Jahre"], b_base$estimates$Garantie["AMCE","Garantie8Jahre"], b_base$estimates$Garantie["AMCE","Garantie10Jahre"])
s_attribute8 <- c(0, b_base$estimates$Garantie["Std. Error","Garantie4Jahre"], b_base$estimates$Garantie["Std. Error","Garantie6Jahre"], b_base$estimates$Garantie["Std. Error","Garantie8Jahre"], b_base$estimates$Garantie["Std. Error","Garantie10Jahre"])

c_attribute9 <- c(0, b_base$estimates$Reichweite["AMCE","Reichweite300km"], b_base$estimates$Reichweite["AMCE","Reichweite500km"], b_base$estimates$Reichweite["AMCE","Reichweite1000km"])
s_attribute9 <- c(0, b_base$estimates$Reichweite["Std. Error","Reichweite300km"], b_base$estimates$Reichweite["Std. Error","Reichweite500km"], b_base$estimates$Reichweite["Std. Error","Reichweite1000km"])

c_attribute10 <- c(0, b_base$estimates$Ladegarantie["AMCE","Ladegarantie1km"], b_base$estimates$Ladegarantie["AMCE","Ladegarantie2km"], b_base$estimates$Ladegarantie["AMCE","Ladegarantie10km"])
s_attribute10 <- c(0, b_base$estimates$Ladegarantie["Std. Error","Ladegarantie1km"], b_base$estimates$Ladegarantie["Std. Error","Ladegarantie2km"], b_base$estimates$Ladegarantie["Std. Error","Ladegarantie10km"])


c <- c(c_attribute1, c_attribute2, c_attribute3, c_attribute4, c_attribute5, c_attribute6, c_attribute7, c_attribute8, c_attribute9, c_attribute10)
s <- c(s_attribute1, s_attribute2, s_attribute3, s_attribute4, s_attribute5, s_attribute6, s_attribute7, s_attribute8, s_attribute9, s_attribute10)

df_b_base <- data.frame(c,s)

b_base$attributes$Tageskarte <- ifelse(b_base$attributes$Tageskarte == "CHF0", "1CHF0", ifelse(
  b_base$attributes$Tageskarte == "CHF150", "1CHF150", b_base$attributes$Tageskarte))
b_base$attributes$Beitrag <- ifelse(b_base$attributes$Beitrag == "CHF0", "2CHF0", ifelse(
  b_base$attributes$Beitrag == "CHF500", "2CHF500", ifelse(
    b_base$attributes$Beitrag == "CHF2000", "2CHF2000", ifelse(
      b_base$attributes$Beitrag == "CHF4000", "2CHF4000", ifelse(
        b_base$attributes$Beitrag == "CHF6000", "2CHF6000", b_base$attributes$Beitrag)))))


df_b_base$att <- c(b_base$attributes$Kaufpreis[1], b_base$attributes$Kaufpreis[2:5], 
                   b_base$attributes$Energiekosten[1], b_base$attributes$Energiekosten[2:5],
                   b_base$attributes$Tageskarte[1], b_base$attributes$Tageskarte[2:4], 
                   b_base$attributes$Verbrennertausch[1], b_base$attributes$Verbrennertausch[2:4],
                   b_base$attributes$Beitrag[1], b_base$attributes$Beitrag[2:5], 
                   b_base$attributes$Ladedauer[1], b_base$attributes$Ladedauer[2:6],
                   b_base$attributes$Unterhaltskosten[1], b_base$attributes$Unterhaltskosten[2:4],
                   b_base$attributes$Garantie[1], b_base$attributes$Garantie[2:5],
                   b_base$attributes$Reichweite[1], b_base$attributes$Reichweite[2:4],
                   b_base$attributes$Ladegarantie[1], b_base$attributes$Ladegarantie[2:4])

df_b_base$policy <- ifelse(df_b_base$att == "CHF20000" |  df_b_base$att == "CHF40000" |  df_b_base$att == "CHF60000" |  df_b_base$att == "CHF80000" |  df_b_base$att == "CHF100000", 1, ifelse(
  df_b_base$att == "CHF2" | df_b_base$att == "CHF5" | df_b_base$att == "CHF10" | df_b_base$att == "CHF20" | df_b_base$att == "CHF30", 2, ifelse(
    df_b_base$att == "1CHF0" | df_b_base$att == "CHF25" |  df_b_base$att == "CHF75" | df_b_base$att == "1CHF150", 3, ifelse(
      df_b_base$att == "CHF0" | df_b_base$att == "CHF50" | df_b_base$att == "CHF100" | df_b_base$att == "CHF150", 4, ifelse(
        df_b_base$att == "2CHF0" | df_b_base$att == "2CHF500" | df_b_base$att == "2CHF2000" | df_b_base$att == "2CHF4000" | df_b_base$att == "2CHF6000", 5, ifelse(
          df_b_base$att == "15min" | df_b_base$att == "30min" | df_b_base$att == "45min" | df_b_base$att == "60min" | df_b_base$att == "75min" | df_b_base$att == "90min", 6, ifelse(
            df_b_base$att == "CHF500" | df_b_base$att == "CHF2000" | df_b_base$att == "CHF4000" | df_b_base$att == "CHF6000", 7, ifelse(
              df_b_base$att == "2Jahre" | df_b_base$att == "4Jahre" | df_b_base$att == "6Jahre" | df_b_base$att == "8Jahre" | df_b_base$att == "10Jahre", 8, ifelse(
                df_b_base$att == "100km" | df_b_base$att == "300km" | df_b_base$att == "500km" | df_b_base$att == "1000km", 9, ifelse(
                  df_b_base$att == "05km" | df_b_base$att == "1km" | df_b_base$att == "2km" | df_b_base$att == "10km", 10, 11))))))))))

df_b_base$policy = factor(as.factor(df_b_base$policy), levels=c(1:11), labels=c("Purchase\nprice","Energy\ncosts","PT day\npass", "Car\nexchange", 
                                                                                "Government\nsubsidy", "Charging\ntime", "Maintenance\ncosts", "Warranty\nperiod",
                                                                                "Driving\nrange", "Charging\nguarantee", "Other"))


df_b_base$att <- factor(df_b_base$att, levels = c(rev(c("CHF20000","CHF40000","CHF60000","CHF80000","CHF100000",
                                                        "CHF2","CHF5","CHF10", "CHF20", "CHF30", 
                                                        "1CHF0","CHF25", "CHF75", "1CHF150", 
                                                        "CHF0","CHF50", "CHF100","CHF150",
                                                        "2CHF0","2CHF500","2CHF2000", "2CHF4000", "2CHF6000", 
                                                        "15min", "30min", "45min", "60min", "75min","90min", 
                                                        "CHF500","CHF2000", "CHF4000", "CHF6000",
                                                        "2Jahre","4Jahre","6Jahre","8Jahre","10Jahre",
                                                        "100km", "300km", "500km", "1000km",
                                                        "05km","1km", "2km", "10km"))),
                        labels = c(rev(c("Baseline: CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000", "CHF 100'000",
                                         "Baseline: CHF 2",  "CHF 5", "CHF 10","CHF 20", "CHF 30",
                                         "Baseline: CHF 0",  "CHF 25",  "CHF 75", "CHF 150",
                                         " Baseline: CHF 0", "CHF 50", "CHF 100", "CHF 150",
                                         "  Baseline: CHF 0", "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                         "Baseline: 15 min", "30 min", "45 min", "60 min", "75 min", "90 min",
                                         "Baseline: CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                         "Baseline: 2 years",  "4 years",  "6 years",  "8 years", "10 years",
                                         "Baseline: 100 km", "300 km",  "500 km", "1'000 km",
                                         "Baseline: 0,5 km", "1 km", "2 km", "10 km"))))




figbuy <- ggplot(df_b_base, aes(x=att, y = c)) +
  geom_pointrange(aes(min = c - 1.95 * s, max = c + 1.95 * s)) +
  theme_bw() + 
  theme(panel.grid.major.x = element_line(size=.3, color="gray80"),
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(size=.3, color="gray80" ),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom",
        strip.background = element_rect(fill = 'white')) +
  coord_flip() +
  ylab("Estimated AMCE") +
  xlab("Attributes") +
  geom_hline(yintercept = 0, lty="dashed") +
  facet_grid(policy ~ ., scales = "free") +
  ggtitle("Conjoint purchase decision")
figbuy

ggsave(figbuy, filename = "./FigureBuy.png", width = 10.5, height = 7.5, unit="in", dpi = 600)


############### Joint graph #############
jGraph <- rbind(df_c_allC, df_b_base)
jGraph$grp <- rep(c("Choice", "Purchase decision"), each=nrow(df_c_allC))

jGraphFig <- ggplot(jGraph, aes(x=att, y = c)) +
  geom_pointrange(aes(min = c - 1.95 * s, max = c + 1.95 * s)) +
  theme_bw() + 
  theme(panel.grid.major.x = element_line(size=.3, color="gray80"),
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(size=.3, color="gray80" ),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom",
        strip.background = element_rect(fill = 'white')) +
  coord_flip() +
  ylab("Estimated AMCE") +
  xlab("Attributes") +
  geom_hline(yintercept = 0, lty="dashed") + 
  facet_grid(policy ~ grp, scales = "free")
jGraphFig

ggsave(jGraphFig, filename = "./FigureJoint.png", width = 7, height = 9, unit="in", dpi = 600)


################## Marginal Means ###############

library(cregg)

df_conj$Tageskarte <- ifelse(df_conj$Tageskarte == "CHF 0", "1CHF 0", ifelse(
  df_conj$Tageskarte == "CHF 25", "CHF 25", ifelse(
    df_conj$Tageskarte == "CHF 75", "CHF 75", ifelse(
    df_conj$Tageskarte == "CHF 150", "1CHF 150", df_conj$Tageskarte))))
df_conj$Tageskarte <- as.factor(df_conj$Tageskarte)
df_conj$Beitrag <- ifelse(df_conj$Beitrag == "CHF 0", "2CHF 0", ifelse(
  df_conj$Beitrag == "CHF 500", "2CHF 500", ifelse(
    df_conj$Beitrag == "CHF 2'000", "2CHF 2'000", ifelse(
      df_conj$Beitrag == "CHF 4'000", "2CHF 4'000", ifelse(
        df_conj$Beitrag == "CHF 6'000", "2CHF 6'000", df_conj$Beitrag)))))
df_conj$Beitrag <- as.factor(df_conj$Beitrag)


#choice
cmm <- mm(df_conj, choice ~ Kaufpreis +  Energiekosten +  Tageskarte +  Verbrennertausch + Beitrag +  
            Ladedauer + Unterhaltskosten + Garantie + Reichweite + Ladegarantie, id = ~ResponseId, design=design, baselines = baselines_base)

plot(cmm)

cmm$policy <- ifelse(cmm$level == "CHF 20'000" |  cmm$level == "CHF 40'000" |  cmm$level == "CHF 60'000" |  cmm$level == "CHF 80'000" |  cmm$level == "CHF 100'000", 1, ifelse(
  cmm$level == "CHF 2" | cmm$level == "CHF 5" | cmm$level == "CHF 10" | cmm$level == "CHF 20" | cmm$level == "CHF 30", 2, ifelse(
    cmm$level == "1CHF 0" | cmm$level == "CHF 25" |  cmm$level == "CHF 75" | cmm$level == "1CHF 150", 3, ifelse(
      cmm$level == "CHF 0" | cmm$level == "CHF 50" | cmm$level == "CHF 100" | cmm$level == "CHF 150", 4, ifelse(
        cmm$level == "2CHF 0" | cmm$level == "2CHF 500" | cmm$level == "2CHF 2'000" | cmm$level == "2CHF 4'000" | cmm$level == "2CHF 6'000", 5, ifelse(
          cmm$level == "15 min" | cmm$level == "30 min" | cmm$level == "45 min" | cmm$level == "60 min" | cmm$level == "75 min" | cmm$level == "90 min", 6, ifelse(
            cmm$level == "CHF 500" | cmm$level == "CHF 2'000" | cmm$level == "CHF 4'000" | cmm$level == "CHF 6'000", 7, ifelse(
              cmm$level == "2 Jahre" | cmm$level == "4 Jahre" | cmm$level == "6 Jahre" | cmm$level == "8 Jahre" | cmm$level == "10 Jahre", 8, ifelse(
                cmm$level == "100 km" | cmm$level == "300 km" | cmm$level == "500 km" | cmm$level == "1000 km", 9, ifelse(
                  cmm$level == "0,5 km" | cmm$level == "1 km" | cmm$level == "2 km" | cmm$level == "10 km", 10, 11))))))))))


cmm$policy = factor(as.factor(cmm$policy), levels=c(1:11), labels=c("Purchase\nprice","Energy\ncosts","PT day\npass", "Car\nexchange", 
                                                                    "Government\nsubsidy", "Charging\ntime", "Maintenance\ncosts", "Warranty\nperiod",
                                                                    "Driving\nrange", "Charging\nguarantee", "Other"))

cmm$att <- factor(cmm$level, levels = c(rev(c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000",  "CHF 100'000",
                                              "CHF 2",       "CHF 5",       "CHF 10",      "CHF 20",      "CHF 30",
                                              "1CHF 0",      "1CHF 150",    "CHF 25",      "CHF 75",     
                                              "CHF 0",       "CHF 50",      "CHF 100",     "CHF 150",
                                              "2CHF 0", "2CHF 500", "2CHF 2'000",  "2CHF 4'000",      "2CHF 6'000",
                                              "15 min",      "30 min",      "45 min",      "60 min",      "75 min",
                                              "90 min",      "CHF 500",     "CHF 2'000",   "CHF 4'000",   "CHF 6'000",
                                              "2 Jahre",     "4 Jahre",     "6 Jahre",     "8 Jahre",     "10 Jahre",
                                              "100 km",      "300 km",      "500 km",      "1000 km",
                                              "0,5 km",      "1 km",        "2 km",        "10 km"))),
                  labels = c(rev(c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000", "CHF 100'000",
                                   "CHF 2",  "CHF 5", "CHF 10","CHF 20", "CHF 30",
                                   "CHF 0",  "CHF 25",  "CHF 75", "CHF 150",
                                   " CHF 0", "CHF 50", "CHF 100", "CHF 150",
                                   "  CHF 0", "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                   "15 min", "30 min", "45 min", "60 min", "75 min", "90 min",
                                   "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                   "2 years",  "4 years",  "6 years",  "8 years", "10 years",
                                   "100 km", "300 km",  "500 km", "1'000 km",
                                   "0,5 km", "1 km", "2 km", "10 km"))))




p_cmm <- ggplot(cmm, aes(x=att, y = estimate)) +
  geom_pointrange(aes(min = estimate - 1.95 * std.error, max = estimate + 1.95 * std.error)) +
  theme_bw() + 
  theme(panel.grid.major.x = element_line(size=.3, color="gray80"),
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(size=.3, color="gray80" ),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom",
        strip.background = element_rect(fill = 'white')) +
  coord_flip() +
  ylab("Marginal Means") +
  xlab("Attributes") +
  geom_hline(yintercept = 0.5, lty="dashed") +
  facet_grid(policy ~ ., scales = "free")
p_cmm

ggsave(p_cmm, filename = "./FigureMMChoice.png", width = 10.5, height = 7.5, unit="in", dpi = 600)



#buy
bmm <- mm(df_conj, buy ~ Kaufpreis +  Energiekosten +  Tageskarte +  Verbrennertausch + Beitrag +  
            Ladedauer + Unterhaltskosten + Garantie + Reichweite + Ladegarantie, id = ~ResponseId, design=design, baselines = baselines_base)

plot(bmm)

bmm$policy <- ifelse(bmm$level == "CHF 20'000" |  bmm$level == "CHF 40'000" |  bmm$level == "CHF 60'000" |  bmm$level == "CHF 80'000" |  bmm$level == "CHF 100'000", 1, ifelse(
  bmm$level == "CHF 2" | bmm$level == "CHF 5" | bmm$level == "CHF 10" | bmm$level == "CHF 20" | bmm$level == "CHF 30", 2, ifelse(
    bmm$level == "1CHF 0" | bmm$level == "CHF 25" |  bmm$level == "CHF 75" | bmm$level == "1CHF 150", 3, ifelse(
      bmm$level == "CHF 0" | bmm$level == "CHF 50" | bmm$level == "CHF 100" | bmm$level == "CHF 150", 4, ifelse(
        bmm$level == "2CHF 0" | bmm$level == "2CHF 500" | bmm$level == "2CHF 2'000" | bmm$level == "2CHF 4'000" | bmm$level == "2CHF 6'000", 5, ifelse(
          bmm$level == "15 min" | bmm$level == "30 min" | bmm$level == "45 min" | bmm$level == "60 min" | bmm$level == "75 min" | bmm$level == "90 min", 6, ifelse(
            bmm$level == "CHF 500" | bmm$level == "CHF 2'000" | bmm$level == "CHF 4'000" | bmm$level == "CHF 6'000", 7, ifelse(
              bmm$level == "2 Jahre" | bmm$level == "4 Jahre" | bmm$level == "6 Jahre" | bmm$level == "8 Jahre" | bmm$level == "10 Jahre", 8, ifelse(
                bmm$level == "100 km" | bmm$level == "300 km" | bmm$level == "500 km" | bmm$level == "1000 km", 9, ifelse(
                  bmm$level == "0,5 km" | bmm$level == "1 km" | bmm$level == "2 km" | bmm$level == "10 km", 10, 11))))))))))


bmm$policy = factor(as.factor(bmm$policy), levels=c(1:11), labels=c("Purchase\nprice","Energy\ncosts","PT day\npass", "Car\nexchange", 
                                                                    "Government\nsubsidy", "Charging\ntime", "Maintenance\ncosts", "Warranty\nperiod",
                                                                    "Driving\nrange", "Charging\nguarantee", "Other"))

bmm$att <- factor(bmm$level, levels = c(rev(c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000",  "CHF 100'000",
                                              "CHF 2",       "CHF 5",       "CHF 10",      "CHF 20",      "CHF 30",
                                              "1CHF 0",      "1CHF 150",    "CHF 25",      "CHF 75",     
                                              "CHF 0",       "CHF 50",      "CHF 100",     "CHF 150",
                                              "2CHF 0", "2CHF 500", "2CHF 2'000",  "2CHF 4'000",      "2CHF 6'000",
                                              "15 min",      "30 min",      "45 min",      "60 min",      "75 min",
                                              "90 min",      "CHF 500",     "CHF 2'000",   "CHF 4'000",   "CHF 6'000",
                                              "2 Jahre",     "4 Jahre",     "6 Jahre",     "8 Jahre",     "10 Jahre",
                                              "100 km",      "300 km",      "500 km",      "1000 km",
                                              "0,5 km",      "1 km",        "2 km",        "10 km"))),
                  labels = c(rev(c("CHF 20'000",  "CHF 40'000",  "CHF 60'000",  "CHF 80'000", "CHF 100'000",
                                   "CHF 2",  "CHF 5", "CHF 10","CHF 20", "CHF 30",
                                   "CHF 0",  "CHF 25",  "CHF 75", "CHF 150",
                                   " CHF 0", "CHF 50", "CHF 100", "CHF 150",
                                   "  CHF 0", "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                   "15 min", "30 min", "45 min", "60 min", "75 min", "90 min",
                                   "CHF 500", "CHF 2'000", "CHF 4'000", "CHF 6'000",
                                   "2 years",  "4 years",  "6 years",  "8 years", "10 years",
                                   "100 km", "300 km",  "500 km", "1'000 km",
                                   "0,5 km", "1 km", "2 km", "10 km"))))




p_bmm <- ggplot(bmm, aes(x=att, y = estimate)) +
  geom_pointrange(aes(min = estimate - 1.95 * std.error, max = estimate + 1.95 * std.error)) +
  theme_bw() + 
  theme(panel.grid.major.x = element_line(size=.3, color="gray80"),
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(size=.3, color="gray80" ),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom",
        strip.background = element_rect(fill = 'white')) +
  coord_flip() +
  ylab("Marginal Means") +
  xlab("Attributes") +
  geom_hline(yintercept = 0.5, lty="dashed") +
  facet_grid(policy ~ ., scales = "free")
p_bmm

ggsave(p_bmm, filename = "./FigureMMBuy.png", width = 10.5, height = 7.5, unit="in", dpi = 600)


############### Joint graph Marginal Means #############
jGraphMM <- rbind(cmm, bmm)
jGraphMM$grp <- rep(c("Choice", "Purchase decision"), each=nrow(cmm))

jGraphMMFig <- ggplot(jGraphMM, aes(x=att, y = estimate)) +
  geom_pointrange(aes(min = estimate - 1.95 * std.error, max = estimate + 1.95 * std.error), position = position_dodge(width=.5)) +
  theme_bw() + 
  theme(panel.grid.major.x = element_line(size=.3, color="gray80"),
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(size=.3, color="gray80" ),
        panel.grid.minor.y = element_blank(),
        legend.position = "bottom",
        strip.background = element_rect(fill = 'white')) +
  coord_flip() +
  ylab("Marginal Means") +
  xlab("Attributes") +
  facet_grid(policy ~ grp, scales = "free")
jGraphMMFig

ggsave(jGraphMMFig, filename = "./FigureJointMM.png", width = 7, height = 9.5, unit="in", dpi = 600)
